Phase diagram and Gap anisotropy in Iron-Pnictide Superconductors 
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Using the fluctuation-exchange approximation, we study an effective five-band Hubbard model 
for iron-pnictide superconductors obtained from the first-principles band structure. We preclude 
deformations of the Fermi surface due to electronic correlations by introducing a static potential, 
which mimics the effect of charge relaxation. Evaluating the Eliashberg equation for various dopings 
and interaction parameters, we find that superconductivity can sustain higher hole than electron 
doping. Analyzing the symmetry of the superconducting order parameter we observe clear differences 
between the hole and electron-doped systems. We discuss the importance of the pnictogen height 
for superconductivity. Finally, we dissect the pairing interaction into various contributions, which 
allows us to clarify the relationship between the superconducting transition temperature and the 
proximity to the ant i- ferromagnetic phase. 



PACS numbers: 

I. INTRODUCTION 



The recent discovery of superconductivity at 26K in 
LaFeAsOi-^F^ opened a new field of a highly inten- 
sive research in the material science^ In a short pe- 
riod of time, the superconducting transition tempera- 
ture T c has been elevated to over 50 K by substitution 
of another rare earths for La, which yields the high- 
est T c outside cuprates. 2 4 At present, there exist the 
1111 systems represented by LaFeAsO, the 122 systems 
with BaFe2As2, 5 the 111 systems with LiFeAs, 6 and the 
11 systems with Fe(Se,Te)i£ These four families hold 
similar Fe-pnictogen layers, and are supposed to pos- 
sess the same superconducting pairing mechanism^ The 
superconducting phase appears in a close proximity to 
the stripe-type anti-ferromagnetic (AF) phase of the un- 
doped systems. Early on, it was argued that the AF 
spin fluctuations originate from the nesting between the 
two-dimensional cylindrical Fermi surfaces (two hole sur- 
faces around T point and two electron surfaces around 
M point), and that they give rise to the sign-reversing 
s- wave (s±) superconducting stated Kuroki et al con- 
structed an effective five-band model Hamiltonian in the 
unfolded Brillouin zone (BZ), which can describe the 
band structure of LaFeAsO near the Fermi level, and an- 
alyzed it within the random-phase approximation (RPA) 
obtaining a similar pairing stated This conclusion was 
confirmed within the RPA studies / 2 ! 13 the third-order 
perturbation theory^ the functional renormalization 



group, 15 and the fluctuation- exchange (FLEX) approx- 
imation j 16 i 17 The s±-wave state has been actively dis- 
cussed as a promising candidate for the pairing symmetry 
in the iron-pnictide superconductors. 

The undoped LaFeAsO parent compound has a stripe- 
type AF ground state. With the electron doping due to 
substitution of F for O, the AF phase abruptly vanishes 
in a first order way, and the superconducting phase ap- 
pears. The transition temperature barely changes with 
further the electron doping over the range x = 0.04 — 0.14 
for LaFeAsOi-^. The NMR-l/Ti study in the nor- 
mal state revealed that the strong AF fluctuations in the 
undoped case are dramatically suppressed with electron 
doping, and a pseudogap behavior was observed above 
x = 0.1 Such a pseudogap behavior has been ob- 
served also with the photoemission spectroscopy^ This 
implies that the pairing mechanism cannot be attributed 
solely to the AF spin fluctuations. The NMR-l/Ti re- 
laxation rate in the superconducting state follows the T 3 
dependence / 8 ! 19 ! 21 but recently, also a T 6 -like behavior 
has been reported. 23 Thus, it is not clear whether the 
power-law behavior reflects the superconducting node or 
not. Rather, the lack of the residual density of states 
suggests a fully-gapped state with a gap minimal 

In BaFe2As2, the electron doping by substitution of 
Co for Fe and the hole doping by substitution of K for 
Ba are available The superconductivity exists in a 
wider region for the hole doping than for the electron 
doping, and even the end material KFe2As2 is super- 
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conducting^ From the NMR-l/Ti study in the normal 
state, the correlation between the AF spin fluctuations 
and T c can be deduced. 27 30 As for the superconducting 
symmetry, in the hole-doped region, T 5 -dependence of 
the NMR-l/Ti (Ref. 13 lh and the exponential behavior 
of the penetration depth 3 - 2 - indicate a fully-gapped super- 
conductivity. This is supported by the direct observa- 
tion with the angle-resolved photoemission spectroscopy 
(ARPES) j 33 i 34 In the inelastic neutron scattering mea- 
surements, development of the resonance peak below T c 
was reported, although whether it means the sign change 
in the superconducting gap or not is not clear yet^~— 
In addition, quite recently, indications were found that 
KFe2As2 is a multi-gap system with line nodeS j 39 i 4Q and 
the As-P system shows a line-nodal behavior with high 

rn 41.42 
± c^~^ 

In FeSe, an enhancement of T c with pressure was re- 
ported, accompanied by an increase in the AF spin fluc- 
tuations ji 3 -! 4 ^ The neutron scattering shows a correla- 
tion between the superconductivity and the stripe-type 
AF spin fluctuations, and below T c , a development of 
a remarkable resonance peak. 4 — The thermal conductiv- 
ity 46 and the scanning tunnel microscope^ 7 - show a fully- 
gapped behavior. Moreover, the phase-sensitive analysis 
seems to be consistent with the s±-wave superconductiv- 
ity. 

Although the s± -state is the prime candidate for ex- 
planation of pnictide superconductivity, whether the ro- 
bustness of the superconductivity with respect to the 
presence of impurities can be understood within the s±- 
scenario is a key issue for the future4£ In addition, the 
recently discovered 42622 system with perovskite-block 
laye r 49 i 50 seems to possess a considerably different band 
structure ^ 51 ! 52 while the T c is comparable with the other 
pnictides. Whether this system possesses the same pair- 
ing mechanism is another open problem. 

In order to understand the phase diagrams and the 
magnitude of T c of iron pnictides, it is necessary to study 
superconductivity and the correlation effects using real- 
istic microscopic Hamiltonians. The thermal-Hall con- 
ductivity 53 and the microwave conductivity 32 in the su- 
perconducting state indicate a strong scattering between 
the quasiparticles, and the mass enhancement factor ob- 
served in the ARPES, 54 the optical spectroscopy* 5 - 5 - and 
the de Haas- van Alphen experiment 56 is as large as 2 — 3. 
In addition, ab initio estimates of the interaction pa- 
rameters suggest that these systems are moderately, not 
weekly, correlated ^2r— 

In this paper, we investigate the superconductivity by 
including the correlation effect within the FLEX approx- 
imation. As in the preceding study (Paper I), 17 we en- 
counter the following problem. In the intermediate corre- 
lation regime, the renormalized band structure deviates 
drastically from the local-density approximation (LDA) 
one, which leads to substantial changes in the Fermi sur- 
face and the magnetic fluctuations, and thus spoils the 
good agreement of the LDA Fermi surface with ARPES 
data. The effect was traced to a shift of the renormalized 



d 3z 2_ r 2 site energy closer to the Fermi level, which leads 
to enlarging of the Fermi surface around the F' in the un- 
folded BZ and shrinking of the other sheets. We believe, 
as discussed below, that it is an unphysical artifact of 
combining the ab initio band structure with the FLEX 
approximation. In Paper I, as a tentative method, we 
shifted the site energy of d 3z 2_ r 2 to preserve the shape 
of the Fermi surface in the renormalized band structure. 
Although this allowed us to investigate the superconduc- 
tivity with the effect of correlations, the value of the shift 
is ambiguous. In Paper II, 60 we constructed and studied 
an effective four-band model for d xy ,d yz / zx and d x 2_ y 2 
assuming that d 3z 2_ r 2 stays below the Fermi level and is 
always irrelevant for the low-energy physics. While the 
Fermi surface is not deformed so much even in the in- 
termediate correlation regime as expected, the magnetic 
structure still drastically changes and the stripe-type AF 
never becomes dominant. This lead us to suggest that 
the high-energy physics such as the interactions between 
localized spins, which is not considered in the FLEX, 
may be important to understand the stripe-type AF in 
iron-pnictides. 

In the present study, we examine the five-band model, 
and propose a simple way to avoid the drastic deforma- 
tion of the Fermi surface by adding a static potential to 
the FLEX self-energy. We argue that such potential mim- 
ics the restoring force due to charge relaxation. With this 
modification we obtain a strong stripe-type AF fluctua- 
tions also in the intermediate correlation regime. Then, 
we evaluate the Eliashberg equation in the intermediate 
correlation regime, and investigate the phase diagram as 
a function of doping. In addition, we clarify how to un- 
derstand the correlation between the AF spin fluctua- 
tions and T c , and where the paring glue for superconduc- 
tivity comes from. 

In the following section, first, we calculate the LDA 
band structure in LaFeAsO, and make ab initio con- 
struction of the effective five-band Hubbard model. In 
Sec MI Al we demonstrate the results of the FLEX for the 
five-band Hubbard model, and discuss what kind of prob- 
lems they have. In Sec lIII Bl we introduce a simple way to 
avoid deformations of the Fermi surface due to electronic 
correlations, which mimics the effect of charge relaxation. 
In Sec lIII Cl we verify that the modified FLEX well works 
even in the inter madiate correlation regime. In Sec lIV) 
with this method we investigate the doping dependence 
of the eigenvalue in the Eliashberg equation. We show 
that the obtained phase diagram qualitatively explains 
the overall feature in LaFeAsO and BaFe2As2, and then, 
the pnictogen height is important for high T c . In ad- 
dition, we investigate the gap anisotropy of the pairing 
function obtained in several doping cases. Finally, we 
suggest that it is oversimplified that the pairing mecha- 
nism in this system is attributed to only the conventional 
AF spin fluctuation. In Appendices A— D, we summarize 
hopping integrals in the five-band model and the techni- 
cal parts in symmetry consideration of Hamiltonian and 
FLEX calculations. 
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FIG. 1: (Color online) The electronic band structure of the 
undoped LaFeAsO. Red lines denote the band structure in 
the LDA calculation. Green dashed lines are the ten band 
model obtained with the MLWFs method. 
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FIG. 2: (Color online) (a) The five-band structure and (b) 
The Fermi surface in the unfolded BZ. Orange and green bars 
in (a) represent weights of 3z 2 — r 2 and x 2 — y 2 orbitals of Fe 
3d, respectively. In (b), the inner square rotated 45 degrees 
denotes the folded original BZ. The shaded triangle stands for 
the irreducible part of the unfolded BZ. 



II. ELECTRONIC BAND STRUCTURE AND 
MODEL HAMILTONIAN 

First, we carry out LDA band-structure calculation us- 
ing WIEN2k package in the APW+local orbital basis 
The band structure of iron-pnictides is known to be sen- 
sitive to the internal coordinate z, which determines dis- 
tance of the pnictogen from the iron layer in the unit of 
the c lattice parameter. Experimentally, the supercon- 
ducting transition temperature T c is well correlated with 
the Fe-X-Fe (X=pnictogen) bond-angle, which is related 
to the internal coordinate z. 62 Thus z is an important pa- 



rameter of the iron-pnictides. 

In our calculations, performed for non-spin-polarized 
undoped LaFeAsO, we use the experimental val- 
ues for the crystallographic parameters [a=4.035 (A), 
c=8.741 (A), z La = 0.1415, and z As = 0.6512]. 1 In Fig.[Q 
we show the calculated band structure. As Ap and O 
2p hybrid bands are located between — 5 ~ — 2eV, the 
Fe 3d orbitals are dominant contribution to the bands 
near the Fermi level, and the La 4/ bands are situated 
around 3 ~ 4eV. The Fermi surface consisting of quasi- 
two dimensional cylinders is essentially the same as the 
one obtained by Mazin et al££ 

Next, we use the recently developed interface^ to 
wannier90 code 65 and construct the maximally localized 
Wannier functions (MLWFs) spanning the Hilbert space 
of the Fe 3d bands. The tight-binding model on the ML- 
WFs basis provides an input for further calculations. The 
tight-binding band structure (see Appendix |A|) , marked 
by the green dashed line in Fig. [TJ represents well the 
LDA Fe 3d bands. Deviations around ~ — 0.5 and ~ 2eV 
originate from strong hybridization with As 4p, O 2p and 
La 5d orbitals. The ten-band model can be unfolded 11 to 
five bands in a doubled BZ, shown in Fig. [2j The band 
structure is similar to that of the preceding studies j n i 17 
The Fermi surface is composed of two hole sheets (a) 
around the r point, two electron sheets (/?) around the 
M point, and a hole sheet (7) around the r' point (the r 
point of the original folded BZ). The a surface and a part 
of the j3 surface are dominated by the d yz and d zx char- 
acter. The 7 surface and another part of the j3 surface 
arise from d x i_ y i band. The Fermi surface is character- 
ized by nesting with Q = (ir, 0), which corresponds to 
stripe-type AF ordering. In fact, the (7r, 0) AF spin fluc- 
tuation is the dominant fluctuation within RPA. 11 This 
fluctuation is predominantly due to scattering between 
the /3 and 7 surfaces with d x 2_ y 2 character. Thus, the 
presence of the 7 surface with the high density of states is 
vital for this AF fluctuation ! 17 i 6Q i 66 Moreover, since the 
(7r, 0) AF fluctuations are considered the principal pair- 
ing glue for the sign-reversing s± pairing, the size of the 
7 surface plays a key role for the high-T c superconduc- 
tivity in this system. However, it has been shown that 
the 7 hole surface is missing in the band structure in the 
iron-pnictides with perovskite-block layer— ^ mentioned 
in Sec. I. 

Construction of the model Hamiltonian is completed 
by adding the on-site Coulomb interaction, 

H ' = § J2 c li* c uv c it* c it° ( ia ) 

it cr 

U' 

+ Y J2 Y C Uv4rn*> C irna>Cila (lb) 
ii^m era' 

^2" ^ ^ C i£a C ima' C ^cr (lc) 
iiy^m era' 

+ y ^ C i£a C i£a C irnaCi m a, (Id) 



where a = ± and a — — cr, and c]^., and are the 
creation and annihilation operators in the basis of real 
harmonics 3z 2 — r 2 , xz, yz, x 2 — y 2 , and xy located at Fe 
sites. Since the Fermi surface of LaFeAsO is quasi- two- 
dimensional cylinder, for simplicity, we restrict ourselves 
to a two-dimensional k z = space. In addition, to meet 
the rotation invariance of the atomic orbitals in the or- 
bital space, we take U = U' + 2 J and J = J' '. 

In the following section, we show the results of the 
FLEX calculations. The technical parts of FLEX are 
summarized in Appendix D. In the actual calculations, 
we take 64 x 64 meshes in the unfolded BZ and 1024 
Matsubara frequencies. In this case, we can safely carry 
out the FLEX calculation for T > 0.002. We set T = 
0.003 throughout this paper. As the numerical analytic 
continuation, we use the Pade approximation. 



III. RESULTS FOR UNDOPED CASE 
A. Band-structure renormalization 

As mentioned in Sec. I, the FLEX calculations for 
the present model encounter severe problems as reported 
in Papers I and I L 17 i 6Q We start with summary and de- 
tailed analysis of the trouble points for the carrier density 
n = 6.00 corresponding to the undoped parent compound 
LaFeAsO. In Fig.[3fa) we show the largest eigenvalue A of 
the Eliashberg equation as a function of U for two choices 
of J = U/6 and U/8. In both cases, the pairing symmetry 
for the maximum eigenvalue is a s±-wave, as in the RPA 
calculations. Although A initially increases with U, for 
J = U/8 it shows a tendency to saturate for U > 
For J = U/8, we cannot reach A = 1 even for large J7, 
that is to say that even at T = 0.003 c± 30K supercon- 
ductivity cannot be realized. While we can obtain A = 1 
for J = U/6 and U > 1.6, the structure of the magnetic 
fluctuations is drastically changed, as discussed in Papers 
I and II. With increasing U the dominant magnetic fluc- 
tuation moves from the stripe-type AF with Qi = (7r, 0) 
into a checkerboard- type AF with Q2 = (ir, 7r), as illus- 
trated in Figs. [3fb) and[3fc). This change is related to 
the renormalization of the band structure. 

In Fig. [3fd), we show the spectral weight p{k,uj) = 
— ^ I m ^iK^ CJ )' obtained from the imaginary part 
of the retarded Green's function, and the quasiparticle 
band £ n fc along the symmetry lines. To evaluate £ n fc the 
self-energy is expanded to the first order in lj, and the 
equation 

Det \(z f k uj - fi)S em - h k lm - ReZ? m (k, 0) | = (2) 



is solved. The mass enhancement factor z 
bital, which is given by 



dZf £ (k,u) 



duo 



for each or- 
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FIG. 3: (Color online) (a) Eigenvalue A as a function of U in 
the FLEX approximation, (b) Spin susceptibility \ s (^O) at 
U = 1.40 and (c) that at U = 1.80 for J = U/6. (d) Spectral 
weight p(k, uj) along the high symmetry lines for U — 1.40 and 
J — U/6. White and orange lines denote the unperturbed and 
the renormalized bands, respectively. 



can be at sufficiently low temperatures replaced with the 
approximate form 



z k 



1 



7rT ' 



(4) 



while ReU£ m (k, iirT) is taken for the real part of the 
retarded self-energv ReZf m (k,0). Thus obtained quasi- 
particle band traces closely the position of the peak in 
the the spectral density [the bright portion in Fig. [3fd)]. 
The quasiparticle bands are strongly renormalized, which 
results in a drastic change in the Fermi surface [see 
Fig. Ufa)]. The Fermi surfaces around r and M shrink, 
and the 7 sheet around F' vanishes completely. 

To understand the effect of this drastic change on the 
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FIG. 4: (Color online) (a) The renormalized Fermi surface for 
U = 1.40 and J = U/6. (b) The orbital-dependent weights 
of the Green's function \Gu(k,i7rT)\ for d 3z 2_ r 2 -orbital, (c) 



z/c^x -orbital, and (d) d x i 2 -orbital. 



dominant magnetic fluctuations, we show the orbital- 
dependent weight of the Green's function \Gu(k, mT) \ in 
Figs. Hfb) -|U(d) . Shown in Fig. HJc), d yz /d zx bands give 
rise to two nesting vectors, Qi = (tt, 0) and Q2 = (?r, 7r), 
as in the non- interacting case. The renormalization of 
the d x 2_ y 2 band leads to relative suppression of the Qi 
fluctuation due to the reduction in the weight around F' '. 
On the other hand, the rise of the d 3z 2_ r 2 band toward 
the Fermi level results in Q2 fluctuations. Altogether, 
the band renormalization leads to strong enhancement 
of the Q2 fluctuation, which then dominates over the Qi 
fluctuation. 

Up to now, however, there is no experimental evidence 
that Q2 fluctuation is strong and that the renormalized 
bands drastically deviate from the LDA bands except 
for effective mass enhancement. Even though, quite re- 
cently, the Fermi surface shrinkage has been reported in 
the vicinity of the AF critical point in As-P system^ 
and thus the tendency observed with the FLEX approx- 
imation seems to be right, 68 the extent of the Fermi sur- 
face renormalization is overestimated. This leads us to 
conclude that the observed changes in the Fermi surface, 
caused by the shifts of the quasiparticle bands, are largely 
artifact of the LDA+FLEX approximation. In partic- 
ular, the deviation connected with the d 3z 2_ r 2 orbital 
is serious. In Paper II, we constructed and studied an 
effective four-band model excluding the d 3z 2_ r 2 orbital. 
However, even in the absence of d 3z 2_ r 2, Qi fluctuation 
is quite suppressed and Q2 channel becomes dominant in 
the intermediate correlation regime. Thus, the straight- 
forward application of FLEX for these models has the 
problem that the system cannot achieve the stripe-type 
AF phase transition. In the following we will argue that 



the Fermi surface should not be substantially changed by 
the correlations and present a simple way to achieve this 
within the present computational scheme. 



B. Modified self-energy and density relaxation 

We start with the empirical observation that LDA is 
surprisingly successful in predicting the Fermi surface ge- 
ometries even in complicated multi-band materials with 
strong electronic correlations such as heavy fermion sys- 
tems jSSiZSl This fact is even more striking when we realize 
that other quantities such as the effective electron mass, 
spin susceptibility, or specific heat may be completely 
wrong. Quite likely, this success of LDA is connected to 
its high accuracy in computing the charge distributions, 
obtained by minimizing the density functional, which 
contains the large electrostatic (Hartree) contribution. In 
real material, the large Hartree term is the main restoring 
force which stabilizes the charge distribution. 

Constructing the effective Hubbard model only the on- 
site interaction within the Fe d shell is treated explicitly 
while the other (large) interaction terms are absorbed 
into the fixed effective site energies and hopping integrals, 
which do not depend on the charge distribution. Thus 
an important feedback mechanism, which stabilizes the 
charge distribution is missing. For example, in the calcu- 
lations of the previous section the renormalized d 3z 2_ r 2 
and d x 2_ y 2 occupancies deviate about 10% from the un- 
perturbed state. This is remarkably large deviation since 
the Fermi surface is small in this system. This leads 
to a remarkable modification of the small sheets of the 
Fermi surface. An obvious solution to the problem of 
the missing feedback is a self-consistent recalculation of 
the effective Hamiltonian for the each FLEX iteration. 
However, it is not feasible with our present computer 
codes. Therefore we take an alternative "poor man's" 
approach. Taking LDA Fermi surface for realistic, we 
restrict its modification due to correlations by subtract- 
ing the static part from the single-particle self-energy. 
Namely, we replace the FLEX self-energy E£ m (k : iuj n ) in 
Ea. (|DTa|) with 

5£ £m (k, iu n ) = Uem(k, iuj n ) - Re£f m (k, 0), (5) 



where Z^ m (fc, iuj n ) is calculated in Eq.(|D4aj) as usual, 
and Re £f m (k,0) is its static energy part obtained by 
analytic continuation to the real axis of them. Although a 
numerical analytic continuation generally includes errors, 
the uo — ^ limit at low temperatures can be obtained with 
high precisions by the Pade approximation using several 
lowest Matsubara frequencies. 



C. FLEX with modified self-energy 

Here, we present the LDA+FLEX results for the un- 
doped model obtained with the modified self-energy. In 
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FIG. 5: (Color online) (a) U dependence of the eigenvalue 
A and the maximum of the spin susceptibility x s (<?,0) in 
the FLEX with modified self-energy, (b) Spin susceptibility 
X s (q,0) at U = 1.20, and (c) that at U = 1.80 for J = U/6. 
(d) Spectral weight uj) along the high symmetry lines for 
U = 1.80 and J = U/6. White and orange lines denote the 
unperturbed band structure rescaled by 1/3 and the renor- 
malized band structure, respectively. 



Fig.[5fa), we show the largest eigenvalue A and the maxi- 
mum value of the spin susceptibility x s (q, 0) as a function 
of U for J = U/6. The pairing symmetry corresponding 
to A remains s±-wave as in IIII A[ Although the value of 
A is almost the same as that in Fig. (3]^a) for U ^1.0, we 
can obtain a monotonic behavior also for the intermedi- 
ate correlation regime, which is different from the erratic 
behavior in Fig. [3^a). (where A is sensitive to a small 
change in parameters.) This is because the electronic 
structure does not change drastically in this method as 
expected. In fact, Fig. [5jb) and [5jc) indicate that the 
structure of the spin fluctuation does not change even for 




FIG. 6: (Color online) U dependence of the mass enhance- 
ment factor for each orbital z\, at (k x ,k y ) — (77r/32,0) and 
(257r/32, 0) with the same parameters as in Fig. Eta). 



U = 1.8, and (7r, 0) fluctuation is robust. Strictly speak- 
ing, the peak position is not commensurate. This draw- 
back of the perfect elimination of the self-energy shift has 
no serious effect in comparison with the drastic change 
from Qi to Q2 spin fluctuations. In Fig. [5jd), we show 
the spectral weight p(k,u) along the symmetry line for 
U = 1.8. The total bandwidth is reduced by about 1/3 
and the renormalization effect is larger near the Fermi 
level, but the Fermi-surface topology does not change. 
In Fig. [6j we present the mass enhancement factors for 
each orbital close to the a and f3 surfaces. The quasipar- 
ticle mass increases toward the AF critical point, which 
is a behavior recently observed in the dHvA experiment 
for the As-P system 56 . 



IV. DOPING AND SUPERCONDUCTING 
PROPERTIES 

In the rest of the paper, we discuss the results obtained 
with the modified FLEX for various dopings of the five- 
band model. 



A. Fermi surface 

We start with the evolution of the Fermi surface with 
carrier doping. Figure [7] shows the Fermi surfaces of the 
five-band model for electron densities in the interval from 
5.52 to 6.16, where n=6.0 corresponds to the undoped 
LaFeAsO. With electron doping, the 7 surface around 
r' point shrinks and vanishes at n ~ 6.12. The disap- 
pearance of this sheet was shown to be important for 
the pseudogap behavior, 17 while its presence is impor- 
tant for the (7r, 0) spin fluctuations ] 17 i 6Q i 66 The electron 
f3 surface around the M point shrinks with hole dop- 
ing. At n = 5.60, it reduces to a pair of Dirac points, 
and with further hole doping it becomes a hole surface 
as in n = 5.52 case. The same evolution is expected in 
(Ba,K)Fe2As2- The Fermi surface at around n = 5.52 is 
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FIG. 8: (Color online) The phase diagram as a function of the 
carrier doping n for U = 1.20 and J = 0.25. Red and green 
lines denote the eigenvalue A for the s±-wave and d x i_ y i- 
wave, respectively. The blue line represents the maximum of 
the spin susceptibility \ s (<??0)- The shading denotes the AF 
phase. 



FIG. 7: (Color online) Evolution of the Fermi surface with the 
carrier doping. Green lines represent the electron surface, and 
the other reddish lines denote the hole surface. The electron 
surface around the M points shrinks into a pair of Dirac points 
at n = 5.60. 



similar to that of the end material KFe2As2^ Thus, the 
systematic calculation for doping dependence is relevant 
for the overall features of the phase diagrams and the gap 
symmetry in the related materials. 



B. Phase diagram 

Now, let us investigate the doping dependence of the 
largest eigenvalue A of the Eliashberg equation. In Fig.[8j 
we show the eigenvalue for s±-wave and d x 2_ y 2-wave 
states, and the maximum value of the spin susceptibil- 
ity X s (<7,0) for U = 1.20 and J = 0.25. There is a re- 
gion of strong AF spin fluctuation on the hole-doped side 
(n < 6.00) of the phase diagram. The s±-wave state 
overall dominates in the proximity of the AF phase. The 
value of A ~ 0.9 suggests that the superconducting phase 
can be reached for strong interaction and/or lower tem- 
perature. The eigenvalue A is remarkably reduced above 
n ~ 6.10 and below 5.60. This corresponds to vanishing 
of the 7 hole surface and the f3 electron surface at these 
points, respectively. The ^2 -wave solution dominates 
for n > 6.16, where the 7 hole surface is absent. The 
s±-wave in this region develops nodes, so-called nodal 
s±-waveSi, as discussed in section ITV CI 

Next, we examine the doping dependence of s±-wave 
eigenvalue for several different parameters. Figure [9^a) 
shows the A vs n doping dependence for various U at fixed 
J = 0.22, and Fig. Efb) for various J at fixed U = 1.20. 
It is intriguing that A is almost flat over a rather wide 
doping range, and sensitive to J rather than U. Increase 
in J enhances the AF spin fluctuation with Q = (7r, 0) 
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FIG. 9: (Color online) Doping dependence of the s±-wave 
eigenvalue (a) for U = 1.40, 1.20, 1.00 at J = 0.22 and (b) 
for J = 0.30, 0.25, 0.22 at U = 1.20. 



and the eigenvalue A. In the large J case, A increases 
when approaching the AF phase boundary, revealing a 
strong correlation between the AF spin fluctuation and 
superconducting T c . Two different behaviors of A can be 
distinguished, (i) for relatively small J, small A insensi- 
tive to carrier doping, and (ii) for relatively large J, large 
A sensitive to carrier doping and the presence of the AF 



TABLE I: Height and internal z coordinate of the As atom. La 
is shorthand of LaFeAsO. Nd and P 100% represents LaFeAsO 
with the pnictogen height corresponding to NdFeAsO and 
LaFePO, respectively. The pnictogen height of P 50% is an 
interpolated value between La and P 100%. 

pnictogen height (A) internal coordinate z 
Nd L38 O6580 
La 1.32 0.6512 
P50% 1.23 0.6408 
P100% 1.14 0.6304 
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FIG. 10: (Color online) Doping dependence of the s±-wave 
eigenvalue for several pnictogen heights with U = 1.20 and 
J = 0.25. The band structures for shorthand Nd, La, P 50%, 
and P 100% have been obtained with the use of the pnictogen 
heights in Table [Q 



phase. These facts are consistent with the doping depen- 
dence of the transition temperature in LaFeAs(0,F) and 
(Ba,K)Fe2As2/Ba(Fe,Co)2As2, respectively^ 

Next, to understand why SmFeAsO and NdFeAsO 
have the highest T c among the iron-pnictides, let us in- 
vestigate the relationship between the pnictogen height 
z and T c , which has been found experimentally, 62 and 
stressed theoretically. 66 We calculate the s±-wave eigen- 
value A for several different values of z, shown in Table HJ 
and the lattice parameters of LaFeAsO, repeating the 
LDA+FLEX procedure described above for the respec- 
tive crystal structures, (see Appendix [XJ) In Fig. [TUJ we 
show the results for U = 1.20 and J = 0.25. We can see 
that A grows with an increasing pnictogen height. In the 
Nd case, A is larger for wider region, especially on the 
electron-doped side, than that in the La case. Figures 
mT a) andfTTlfb) show the total weight of the Green's func- 
tion \Qu{k^mT)\ and the spin susceptibility x s (q^) 
at n = 6.08. The bright parts in the former represent the 
Fermi surface. With the decreasing pnictogen height, 
the 7 hole surface around F' point visibly shrinks and 
becomes dim. Correspondingly, the (tt, 0) spin fluctua- 
tion is suppressed. Thus, the pnictogen height is very 
important for the high T c since it controls the size of the 
7 hole surface and the magnitude of the spin fluctuation. 
The observed trend is consistent with the experimental 
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FIG. 11: (Color online) (a) Contour plot of the total weight 
of the Green's function Yt \Su(k, i7rF)\ and (b) Spin suscep- 
tibility x s (<3f)0) at n = 6.08 for several pnictogen heights. 

data showing that NdFeAsO and SmFeAsO hold high T c 
even for heavily electron doping. 



C. Gap anisotropy 

In this section, we examine anisotropy of the gap func- 
tions for various dopings. To this end, we plot the band 
diagonal anomalous Green's function, 

T nn (k,i7TF) = (^ n )X(fc,^)i, (6) 

ira 



obtained from Eq. (|Dlc|) with the eigenvectors of the 
Eliashberg Eq. (|D15|) . It singles out the gap amplitude 
on the Fermi surface. 

In the hole-doped region, the strong (jr, 0) spin fluc- 
tuations render the s±-wave a likely gap function. For 
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FIG. 12: (Color online) Spin susceptibility \ s 0) ( a ) a ^ 
n = 5.76 and (b) at n = 6.12 for U = 1.20 and J = 0.25. 



instance, at n = 5.76, the (7r, 0) spin fluctuation is re- 
markably enhanced as shown in Fig. [T2l a), and the gap 
function has no nodes as shown in Figs. [T3t a) and fl3T b). 
In this case, the gap function on the a surface has almost 
the same amplitude as that on the j3 surface (with op- 
posite sign), and about a half of the 7 surface one. This 
is consistent with the gap structure observed in ARPES 
on (Ba,K)Fe2As2^^ The ratio of the a to 7 surface 
gap amplitude changes gradually with carrier doping, as 
illustrate by n = 5.60 in Fig. [T4Ta). 

In the electron-doped region, the AF spin fluctuations 
are suppressed as shown in Fig.[12fb). While the s±-wave 
pairing symmetry is still favorable, the gap function be- 
comes remarkably anisotropic on the /3 surface, as shown 
in Fig. fT3f d). Further electron doping leads to a sign 
reversal, i.e., it becomes so-called nodal s±-wave shown 
in Fig. fT3Tf). However, the corresponding eigenvalue is 
small, and the d x 2_ y 2-wdive pairing prevails. Thus we 
may expect a fully-gapped state in the hole-doped re- 
gion, and gap minima or line-nodes in the electron-doped 
region. This may be the key to understanding of the 
material-dependent nodal features as the fully-gapped 
behavior in (Ba,K)Fe2As2, and the nodal/nodeless be- 
haviors in Ba(Fe,Co)2As2 and LaFeAs(0,F). 

Let us take a detailed look at the gap anisotropy. Even 
at n = 5.76, the gap amplitude around the M point 
exhibits some anisotropy. Comparing the gap amplitudes 
at the points marked a and b in Fig. fT3T b), we find that 
the amplitude at a is the same as that on the a surface, 
while the amplitude at b is about the same as on the 7 
surface. This reflects the fact that the Fermi surface at 
a is dominated by the d yz /d zx orbitals, and that at b 
by the d x i_ y i orbital. Thus, the gap amplitude on the 
f3 surface will be always anisotropic. The appearance of 
the distinct gap minima on the j3 surface in the electron 
doped region is related to the shrinkage of the 7 surface. 

Finally, let us comment on the end material of the 122 
series, KFe2As2- Quit recently, it has been reported that 
it exhibits a two-gap nodal behavior. 3 -^ In our model, 
this material corresponds to a hole doping of 0.5, close 
to the filling n = 5.60, which leads to a Dirac-cone band 
structure. Since the gap functions for these two fillings 
barely differ, let us discuss only the n = 5.60 case, in 
which the unrenormalized f3 surface reduces to a point. In 




M M 

FIG. 13: (Color online) Band-diagonal anomalous Green's 
function J nn (fc,i7rT) for the s±-wave at n = 5.76, 6.12 and 
6.20. (a), (c), (e) correspond to the anomalous Green's func- 
tion on the third band, and (b), (d), (f) on the second band^. 




M M 

FIG. 14: (Color online) Band-diagonal anomalous Green's 
function J r nn (k, inT) at n = 5.60. (a) and (b) correspond to 
the third and the second bands in the s±-wave case, respec- 
tively, (c) and (d) in the d x 2_ y 2-w&ve case. 
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the strong coupling theory of superconductivity, electrons 
within a finite energy window around the Fermi level 
participate in the pairing, and the anomalous Green's 
function can have a large amplitude even at k points far 
from the Fermi surface. Indeed, we find a large weight 
in a finite region around the M point, see Figs. [T4l b) 
and [T4lfd). Another interesting point is the small mag- 
nitude of the gap function on the 7 surface relative to 
that around the r and M points, both in the s±-wave 
and 6^2 -wave states. Observation of two gap of dis- 
tinct sizes is consistent with the recent experiments in 
KFe2As2^ Nevertheless, this conclusion must be taken 
with care. In addition to the small value of A, we overesti- 
mate the <i 3z 2_ r 2 -orbital contribution in the heavily hole- 
doped region. In the present two-dimensional model, the 
c-axis dispersion due to the d 3z 2_ r 2 orbital is neglected. 
Although it is not crucial as long as the d 3z 2_ r 2 band is 
far from the Fermi level, its presence close to the Fermi 
level suppresses the eigenvalue of the d £C 2_ 2/ 2-wave rela- 
tive to the s±-wave. Therefore, we expect the d x 2_ y 2- 
wave eigenvalue to be comparable with the s±-wave one 
in a more realistic three-dimensional calculation. 



D. Pairing Glue 




In this section, we analyze the paring mechanism. Al- 
though we have obtained the phase diagram insensitive 
to carrier doping and the presence of the AF phase in rel- 
atively small J case, we need to clarify what is the glue 
for the superconducting pairs in this case. In Fig. [15j we 
investigate the J = 0.22 and 0.25 cases. Within the RPA 
and FLEX formalism the pairing interaction can be sep- 
arated into contributions from the spin and charge sec- 
tors. The former includes all spin fluctuations described 
by Xi2,34> an d the latter consists of all charge/orbital 
fluctuations contained in X12 34- ^ n Fig. [15l the "spin 
sector" denotes contribution of X12 34 only? i-e., we omit 
X12 34 when solving the Eliashberg equation. Clearly, the 
contribution of the charge sector is negligible. 

Next, we consider only the contribution from the 
orbital-diagonal spin fluctuations Xu mm- These add 
up to the spin susceptibility Xdiag = ((^i ' 1 = 
J2i m Xii mm' wmcn describes the correlations between 

spin dipole moments Sf = Y^i, a p c \ia G ap c wP- In Fi g-HH1 
we denote this contribution with "spin dipole." We find 
that for both choices of J, the contribution of the spin 
dipole-dipole fluctuations reflects the proximity of the 
AF phase, that is, the eigenvalue A overall increases 
towards the AF phase boundary. Such a behavior is 
rather remarkable in the light of flatness of the full 
A vs n. Moreover, the dipole-only eigenvalue is only 
60 — 70% of the total value of A. The rest comes from 
the orbital-off-diagonal spin fluctuations, which represent 
correlations between higher-order multipoles, especially, 
Xu,£m + Xti y m£ as suggested in Paper I, This gives a cou- 
pling between spin-dipole and spin-quadrupole. Thus, 



FIG. 15: (Color online) Phase diagram for (a) J = 0.22 and 
(b) J — 0.25 at U — 1.20. The violet line represents the max- 
imum of the spin susceptibility X s (q^)- The red, orange, 
and green lines denote the eigenvalue A with the total con- 
tribution, the spin sector and just the conventional AF spin 
fluctuation, respectively. 

although the conventional AF spin fluctuations provide 
the largest contribution to the paring glue, the higher- 
order multipolar fluctuations given by the off-diagonal 
elements of Xi2,34 assist the superconducting pairing, and 
push the transition temperature up. Therefore, the cor- 
relation between T c and the proximity of the AF phase 
is seemingly weak. 



V. SUMMARY AND CONCLUSIONS 

Using a combination of ab initio band structure and the 
FLEX approximation, we have constructed a five-band 
model of iron pnictides and studied the doping depen- 
dence of superconductivity. We have employed a simple 
procedure to fix the shape of the Fermi surface to its LDA 
shape. 

We have found that the superconductivity is stable in 
a wider interval for hole doping than for electron dop- 
ing. The s±-wave state corresponds to the largest eigen- 
value of the Eliashberg equation A over a wide doping 
range, only in the heavily electron doped region d-wave 
state becomes more favorable. For a relatively small 
Hund's coupling J, the eigenvalue A is relatively small 
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and insensitive to carrier doping. On the other hand, 
for large J, A is large and more sensitive to the carrier 
doping and the proximity of AF phase. These obser- 
vations are consistent with the behavior of the super- 
conducting transition temperature T c in LaFeAs(0,F) 
and (Ba,K)Fe2As2/Ba(Fe,Co)2As2, respectively. Fur- 
thermore, to understand why the highest T c is found in 
SmFeAsO and NdFeAsO, we have investigated the influ- 
ence of the pnictogen height on A. The pnictogen height 
has a great impact on the size of the Fermi surface around 
the r' point formed by the d x i_ y i orbital. In the Nd 
case, the F' sheet is large resulting in a large A over 
a wider doping region, especially on the electron-doped 
side, than in the La system. This agrees with the fact 
that SmFeAsO and NdFeAsO hold the high T c even when 
heavily electron doped. 

As for the anisotropy of the gap function, we find 
a fully-gapped state in the hole-doped region, and re- 
markably anisotropic gap function around the M point 
in the electron-doped region. This may be the key to 
understanding of the material-dependent nodal behav- 
ior, such as fully-gapped behavior in (Ba,K)Fe2As2, and 
nodal/nodeless behavior in LaFeAs(0,F). Furthermore, 
we have found indications of gaps of two distinct sizes 
for the end 122 material KFe2As2- 

Finally, concerning the pairing mechanism, we have 
explained why the correlation between T c and the pres- 
ence of the AF phase is seemingly weak in this system. 
Only 60 — 70% of the total pairing interaction origi- 
nates from the diagonal components of spin fluctuation, 
which corresponds to the conventional spin-spin corre- 
lation. The remaining part originates from correlations 
involving higher-order spin multipoles. This additional 
pairing glue naturally comes from the multi-band char- 
acter of iron-pnictides with several different orbital con- 
tributions at the Fermi surface. Therefore, it is an over- 
simplification to attribute pairing in these materials only 
to the conventional AF spin fluctuations. 



Appendix A: Hopping integrals 

Let us here provide a set of in-plane hopping integrals 
on the MLWFs basis, te m [Ax, Ay], for each band calcu- 
lation in Table HI where [zAx,/A^] denotes the in-plane 
hopping vector, and (Jim) the orbit als. Note that the 
x- and ^-axes point toward neighboring Fe atoms while 
x and y in the orbitals (1 : d 3z 2_ r 2, 2 : d xz , 3 : d yz , 
4 : d x 2_ y 2, 5 : d xy ) are those in the coordinate system for 
the original unit cell. In Table HH we list t£ m [Ax, Ay] up 
to the fifth neighbors, the magnitude of which is larger 
than 0.005eV. We can obtain the principal hopping in- 
tegrals for each band structure from this table and the 
relation t£ m [Ax, Ay] = t m i[-Ax, -Ay]. 



Appendix B: Symmetry Consideration for the 
hopping matrix 

We here give careful consideration to symmetry oper- 
ations in the two-dimensional model Hamiltonian. The 
unperturbed Hamiltonian is defined as follows: 

H = Y1 h £mCi £a Ckma (Bla) 
k£ma 

= ^mAla^kma (Bib) 

k£ma 

= ^2e kn al na a kn(T , (Blc) 

kna 

where cj^ a , and c;, mcr are the Fourier transforms of c\ icT , 
and Ciraa. The hopping matrix h^ m has the following 
form: 
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where " i " is the imaginary unit, and " e " and " o " 
denote, respectively, even- and odd-parity functions of 
the wave vector k. Introducing a sign function with s = 
+/— below/above the line k x = k y in Fig.(2fb), io = zs|o| 
with | o | = so, which is unchanged under the inversion. 
In this case, using the transformation 
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TABLE II: Hopping integrals Um[Ax,Ay] for Nd, La, P 50%, and P 100% in Table [J in units of eV. a y , I, and a d mean 
tem[Ax,—Ay], tim[— Ax, — Ay], and tim\Ay, Ax], respectively, "zb" and "±(£ f m f ) n in the columns denote ±.ttm[Ax, Ay] and 
=bt^/ m / [Ax, Ay], respectively. 
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can be cast into a real symmetric form, 
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which is unchanged under the inversion. The third line 
of Eq. (|Bl[) shows the diagonal form obtained from h^ m 
by the unitary transformation u^ n . The band indices n 
are set according to the main orbital component, not in 
the order of energy eigenvalues. In this case, the unitary 
matrix u 1 ^ nn has the same matrix form as the above h^ m , 
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With the use of this unitary matrix, the direct transfor- 
mation from Ckmcr to dkna can be denned by 



' ( u £n) a kna ~ ^ki^in^kn^ 



"Mo- 



Ill this case, the relation 



kno 



(^mn)* = U 



k 

mm 



(B7a) 
(B7b) 

(B8) 



corresponding to the time-reversal symmetry of the 
Hamiltonian naturally holds. Thus, we can fix the com- 
plicated phase factors accompanying the diagonalization 
at each k point. In addition, considering the mirror sym- 
metry about the two lines k x = k y and k x + k y = 2tt 
through the F' point, we can carry out any numerical 
calculations in the reduced zone as shown in Fig. [2fb), 
which improves the accuracy and speeds up the calcula- 
tion. 
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Appendix C: Green's functions 

We here summarize several generic relations for the 
Green's functions. First, let us define the normal and 
anomalous Green's functions for orbitals i and m. We as- 
sume a paramagnetic normal state, and the spin-singlet 
symmetry for the superconducting state. Then the 
Green's functions at the wave vector k and the imagi- 
nary time r are given by 

& m (fc,r) = ((^WcL(0))), (Cla) 
oT grn (k,T) = ((c fc ^(r)c_ fcm ^(0))), (Clb) 

vflmfar) = ((cl^(r)cL.(0))). (Clc) 

Here, ((A(r)B(O)}} = - (T T [A(r) B (0)]} with the conven- 
tional notation. From the above definition and the time- 
reversal invar iance, we obtain the following relations: 

Gme(k,r)* = Gem(k,r) =G£m(-k,r)*, (C2a) 
Fmi(-k,-T)= Fimik^r) =Fim(-k,Ty, (C2b) 

= J^(-fc, r) = 7j m (fc, -r) = jj m (-fe, -r)*. (C2c) 

By the Fourier transformation from k to r, these relations 
are rewritten as follows, 

Gme(-r,r)* = Gim{r,r) =Gem(r,T)*, (C3a) 
^Vt^C-^j-t) = Tim{r,r) =J r £m(r : ry, (C3b) 

= 4^"^) =^L(^-^)=4m(^-^r- (C3C) 

Thus, Gtmi?, T ) an d r ) are rea l functions. On the 

other hand, by the Fourier transformation from r to the 
fermionic Matsubara frequency uo n = (2n + l)7rT, we can 
obtain 

Gmi(k, -iw n )* = Gim(k, iu n ) = Gim(-k, -iu n )* , (C4a) 
Fmi(-k, -iu n ) Fim(k, iu n ) Ji m (-fc, -icj n )*, (C4b) 

= Fml(~ k i iuj n) = FirJki-iUn) = J"j m (-fc, W*>n)*- (C4c) 

Concerning the wave vector fc, each component of the 
normal Green's functions holds the same irreducible rep- 
resentation in the space group as the corresponding hop- 
ping matrix h^ m . Depending on the parity even or odd, 
the last equality in Eq. (jC4aj) becomes 

Gim(k, iuj n ) = ±Gim(k, -iu n )* . (C5) 

By Eqs. (|B3|) and (|B4|h it is convenient to introduce 
Green's functions, 

Gim(k,iuJn) = / dre lUJnT {{c Ma {r)c ] krn(7 )), (C6a) 
Jo 

f P 

J r £m(k,lU n )= / dre %UJnT '((Ckia^C-kma)), (C6b) 

Jo 



in Cfe mcr representation with the inverse temperature /? = 
1/T. In this case, 

(k,iuj n ) i£ m , (C7a) 
Ji m (fe, icj n ) = Tnm(k, iuj n ) i fc * m , (C7b) 
with i_ fcm = i fc * m , and then, 

& m (fe, zu; n ) = & m (fc, -zcj n )* = ^(fe, zo; n ), (C8a) 

ftm(k, iu n ) = ^£m(k, -iu n )* = T m £(k, iu n )* . (C8b) 

Thus, we can eliminate the sign depending on the par- 
ity. In addition to these relations, a value at k point 
in Green's function holds a simple relation to that at the 
star of k under the space-group symmetry. Thus, Green's 
functions at all k points in the unfolded BZ can be gen- 
erated from those values at k points in the reduced zone. 
These relations obtained in this section are practical in 
actual calculations. 

Appendix D: FLEX formalism 

We here summarize the formulation of FLEX, following 
Ref. 73 . The linearized Dyson- Gorkov equations for spin- 
singlet pairing with the abbreviation k = (k,iu) n ) are 
given by 

Gim(k) = GL(k) + G° it (k)U iw (k)Gm'm(k), (Dla) 

Flrn{k) = G££> (k)Gmm> (-k)A£>m> (k) (Dlb) 

= Gee f (k)Gmm'(k)*Ae'm'(k), (Die) 

where Ei m (k) and A £m {k) are the normal and the 
anomalous self-energies, respectively. In Cfc mcr represen- 
tation, these equations can be rewritten as 

Glm(k) = G° £m (k) + G° U >{k)Z£>rn>{k)Grn>rn{k), (D2a) 

^£m(k) = Gee f (k)Gmrn'(k)*£e f rn'(k), (D2b) 

where Z^ m (/c) and A£ m (k) are, respectively, defined by 

^£m(k) = ik £ ££m(k)ikm, (D3a) 
Aim{k) = i k£ A £m (k)i km . (D3b) 

In the FLEX approximation, the self-energies are de- 
scribed with effective interactions composed of ladder 
and bubble diagrams as follows, 

^£m(k) = ^ Vw,mm' (q)Gl>m> (k ~ q), (D4a) 

q 

A lm {k) = -J2 V &',m>m(q)^m<(k - q), (D4b) 

q 

where the effective interactions Vw ,mm' (q) and 
Vt £ , jmm ,(q) are, respectively, given by (££' \mm r ) el- 
ement of matrices, V(q) and V s (q). Here, v n = 2nnT 
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in q = (q,iv n ) is the bosonic Matsubara frequency. 
In the present study, for simplicity, let us consider 
only particle-hole processes and omit particle-particle 
ladder processes in these effective interactions. This 
simplification can be expected to be justified not only 
qualitatively but also semi-quantitatively in the case 
where the spin fluctuation dominates. In this case, the 
effective interactions are, in the matrix form, given by 



lows: 



(D5a) 



V(q) = U n - 2U n - Unxo(q)u n , 
V s (q) = U n + \tj s x s {q)U s ~ \u c X c (q)U c , (D5b) 



with the bare vertices U s,c — =p Uf\- and three sus- 
ceptibilities, Xo(q), X s (q) and X c (q). Each component of 
the bare vertices is defined by 



(U^i)ee,mm = U' 



(D6a) 
(D6b) 
(D6c) 
(D6d) 



E 



vte) = -^2Gem(k + q)Gm>e>(k) 
k 



(DlOa) 



(DlOb) 

The second line can be obtained with the use of the 



Fourier transformation and the relation of Eq. (|C3a|) . In 
these susceptibilities, a simple relation between the upper 
and lower triangular components holds. This considera- 
tion also saves a memory in actual numerical calculations 
and speeds up the calculations. 

In the FLEX approximation, first of all, we evaluate 
the eigenvalue e n u and the unitary matrix u^ m for the 
unperturbed Hamiltonian, h^ m of Eq. (|B5|) , and then, cal- 
culate the chemical potential \i by the condition that the 
electron density is a given n, 



E 

nka 



f (ink) =2^! (ink), (Dll) 



where £ n fc = e n fc — /i, and /(e) = l/(e^ e + l) is the Fermi- 
Dirac distribution function. With the use of the obtained 
unitary matrix £tf m , the noninteracting Green's functions 
are given by 



(U^)u,mm — U' — J 
(U^)em,£m = J — U' 



(D7a) 
(D7b) 



where i ^ m and the other components are zero. 

Susceptibilities X s (q) and x c (q) represent, respectively, 
susceptibilities for spin sector and charge sector, which 
include the RPA-like enhancement of the irreducible sus- 
ceptibility Xo(q). 

X s (q) = MO) + Xo(Q)U s X s (q), (D8a) 
x c (q) = Xo(q) - xo(q)u c x c (q)- (D8b) 

These susceptibilities contain all informations for not 
only spin and charge/orbital fluctuations but also higher- 
order multipolar fluctuations. The specific matrix form 
of Xo(q) is denoted as 74 
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ii /xSi,n 



xo(q) 



21 



12 



A21,ll 



Xl2,ll 



21 

v° 

All, 21 
v° 

A21,21 



Xl2,21 



12 



All, 12 
v° 

A21,12 



Xl2,12 



V 



(D9) 



/ 



With the use of the normal Green's functions, each com- 
ponent of the irreducible susceptibility is defined as fol- 



mn 



1 



(D12a) 



g° em (k,r) = -u k ln u k mn (\ - f(i nk ))e-^. (D12b) 

We next transform Q^ m (k^r) into ^ m (r,r), and then, 
evalua te Eqs.([D8]) via Eq. (|D10bj) . From Eqs. (|D4aj) and 
(jD5a|) , we can obtain the normal self-energy. With 
Eqs. (|D2a|) and (|D3a|h we obtain new normal Green's 
functions, determining a new chemical potential as 

n = 2 ^ f(€nk) + 2 ^ (pmm(k, iuj n ) ~ 5L( fc ? ^n)) • 
nk rnk 

(D13) 

We carry out this self-consistent procedure until relative 
errors in Z}£ m (k) becomes less than 10 -4 . In this paper, 
we define the magnetic transition with x s 0) > 100, 
where X s (q) is the spin susceptibility defined by 



X' 



\q) = J2 x kmm(q)- 



(D14) 



tm 



Concerning the superconducting transition, from 
Eqs.jnic]) and ([DIE]) , we obtain the Eliashberg equa- 
tion, 

k' (D15) 

X^(*')5m"m'(*')^'m'(fc / ) 

with the eigenvalue A. The transition temperature T c 
can be obtained as the temperature when the maximum 
eigenvalue A is unity. 
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